A framework for correcting the noise-induced bias in noisy magnitude MR signals

INTRODUCTION MR signals are complex numbers where the real and imaginary components are independently Gaussian distributed [1]. The phase of the complex MRI signal is highly sensitive to many experimental factors, e.g., see [1,2], and as such, the magnitude of the complex MR signal is used instead in quantitative studies. Although several techniques have been proposed to correct the phase error [2], the magnitude of the complex MR signal (hereafter, magnitudeMR signal) remains the most commonly used measure in MRI. While the magnitude MR signal is not affected by the phase error, it is not an optimal estimate of the underlying signal intensity when the signal-to-noise ratio is low [1] because magnitude MR signals follow a Rician distribution [3] rather than a Gaussian distribution. Although several correction methods have been proposed [1,3-6] to address the noise-induced bias in magnitude data, these methods do not produce corrected data that are Gaussian distributed. Here, we present a scheme to remove the noise-induced bias in noisy magnitude MR signals by making noisy Rician signals Gaussian-distributed. METHODS A simple example illustrates the idea behind the proposed framework: suppose the noisy magnitude signals are drawn from a family of Rician distributions all of which are characterized by different location parameters but with the same scale parameter (e.g., diffusion-weighted signal as a function of b-value or fMRI signal as a function of time). The proposed framework attempts to transform the noisy magnitude signals such that each of the noisy transformed signals may be thought of as if it were drawn from a Gaussian distribution with different mean but the same standard deviation. There are three stages in the proposed scheme. First, a data smoothing method (one, two, or higher-dimensional methods) is used to obtain the average values of the noisy magnitude signals (a penalized spline model [7] is used in this work). The degree of smoothness is selected based on the method of generalized cross-validation [8]. Second, a novel iterative method, different from, but motivated by [6], is employed to take both an estimate of the average value of a noisy magnitude signal and an estimate of the standard deviation of the Gaussian noise, obtained from the image background [1], to an estimate of the average value of the underlying signal intensity. Third, the corresponding noisy Gaussian signal of each of the noisy magnitude signals is found through a composition of the inverse cumulative probability function of a Gaussian random variable and the cumulative probability function of a Rician random variable. The third stage is exactly a Gaussian random number generator if the input data are Rician-distributed. RESULTS We illustrate the performance of our approach on an excised rat hippocampus data set acquired in a 14.1T narrow-bore spectrometer with a pulsed gradient stimulated echo pulse sequence. The imaging parameters were: TE=12.6ms, TR=1000ms, resolution=(78x78x500)μm, matrix size=(64x64x3), number of repetitions=4, diffusion gradient pulse duration (δ)=2ms, and diffusion gradient separation (Δ)=24.54ms. The data set contains a total of 33 images with different diffusion gradient strengths increasing from 0 to 2935mT/m in steps of 91.75mT/m. One diffusion weighted image is shown in Figure 1A. Two neighboring pixels indicated with a red square were selected for further analyses. The noisy magnitude signals of each of the pixels as a function of b-value are shown in Figs. 1B-1C as red dots. The red curve in each of the panels is obtained through a least squares fit of a bi-exponential function to the noisy magnitude signals. We applied the proposed scheme on the noisy magnitude signals (red dots); the resultant or transformed signals are displayed as blue dots in Figs. 1B-1C. The blue curve in each of the panels is obtained through a least square fit of a bi-exponential function to the noisy transformed signals (blue dots) based on the proposed framework. Note that the penalized spline with a truncated polynomial basis of degree 4 and with 4 knots was used in this example. The estimated Gaussian noise standard deviation was 0.88. If both the estimated Gaussian noise SD and each of the blue curves are assumed to be the ground truth values then the expected value (or the first moment) of a Rician distribution as a function of b-values can be computed and is shown in dark gray; these expected values are in good agreement with the red curve, which is an indication that the blue curve is a good approximation of the underlying signal intensities. Note that the increase in variability in the transformed signals at low signal-to-noise ratio is theoretically known and is not unexpected [6]. DISCUSSION & CONCLUSION The proposed scheme is general and is not restricted to diffusion MRI. It can be used in fMRI and other quantitative MRI methods by selecting different smoothing models. For example, the proposed method can be used to correct the noise-induced bias in high angular resolution diffusion MRI or diffusion tensor data by using the spherical harmonic spline model [11] in the first stage of the proposed scheme. The proposed scheme is the first method capable of obtaining corrected data that are distributed evenly in both the positive and negative axes when the signal-to-noise ratio is very close to zero, which is a very important but simple criterion for testing the accuracy or lack thereof of a correction scheme. We should point out that none of the previously published methods [1,3-6] satisfies this criterion because these methods cannot produce corrected data that have negative values. It is important to note that the transformed signals are of interest here rather than the Rician signals because the transformed signals are Gaussian-distributed and are generally more amenable to statistical analysis. As is well known in quantitative MRI, many anatomically or physiologically relevant parameters are usually estimated from models assuming Gaussian-distributed signals. To conclude, the proposed scheme is a practical and effective method for removing the noise-induced bias in noisy magnitude MR signals. The present approach is a major advance in facilitating and improving all subsequent data analysis and processing steps in a quantitative MRI pipeline [9]. REFERENCES [1] Henkelman RM. Med Phys 1985; 12: 232-233. [2] Liu et al. JMR 1990; 86: 593-604. [3] Gudbjartsson et al. MRM 1995; 34: 910-914. [4] McGibney et al. Med Phys 1993; 20: 1077-1078. [5] Miller et al. MRI 1993; 11: 1051-1056. [6] Koay et al. JMR 2006; 179: 317-322. [7] Ruppert et al. Semiparametric regression: CUP; 2003. [8] Wahba G. Spline models for observational data: SIAM; 1990. [9] Zijdenbos et al. IEEE TMI 2002; 21: 1280-1291. Fig. 1. (A) A diffusion-weighted image of a hippocampus with a red square indicating two neighboring pixels selected for further analyses. The data and results are shown in (B) and (C). In each of the figures (B and C) above, the red points are the noisy magnitude signals and the blue points are the corrected signals obtained through the application of the proposed scheme on the red points. Each of the red curves is a smoothed curve obtained through a bi-exponential fitting to the noisy magnitude signals while each of the blue curves is a smoothed curve obtained through a bi-exponential fitting to the transformed noisy signals obtained through the proposed scheme.